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Research Accomplished 


Analysis for Delamination Initiation in Postbuckled Dropped-ply Laminates 
Objective and background 

The objective of this study is analyze the geometrically nonlinear response and delamination 
initiation of compression-loaded dropped-ply laminates. Specimens made from AS4/3502graph- 
ite-epoxy tape were tested in axial compression in an earlier study (Curry, Johnson, and Star- 
nes, 1986) and a schematic of these specimens is shown in Fig. 1. The specimens differed in width 
b and in the layup of the dropped plies. Specimen stacking sequences were 
[ (±45/0/90) s , D, (±45/0/90) s ] in which a zero degree lamina is parallel to the x-axis 

(load axis) and D denotes the layup of the dropped plies. Hence, the layup of the thin section of 
the specimen was fixed to be a quasi-isotropic laminate sixteen plies thick. The thick section, 
which contained the dropped-ply sublaminate, varied from eighteen to twenty-four plies thick. 
Specimens exhibiting a geometrically nonlinear response to delamination initiation had width b 
equal to 3.0 in (7.62 cm). 
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► 
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Fig. 1 Schematic of AS4/3502 graphite/epoxy dropped-ply laminated plates tested in axial 
compression by Curry, et al., 1986. 


Recent publications 

Work continued on two publications that document the analyses for those specimens that 
delaminated at the dropped ply site. The references for, and the status of these papers are: 

• Davila, Carlos G., and Johnson, Eric R., “Analysis for Delamination Initiation in Post- 
buckled Dropped-Ply Laminates,” AIAA Journal, Paper J20239, accepted for publication 
July 25, 1992. 

• Johnson, Eric R., and Ddvila, Carlos G., “Geometrically Nonlinear Transition Elements 
for Laminate Models Containing Interlaminar Stress Risers,” Composites Engineering, 


NASA Grant NAG- 1-537 


2 



submitted July 27, 1992 and under revision, 

Copies of these publications have been sent to the technical monitor. 

Analysis of the [90 2 ] dropped-ply specimen 

The compression test of a [9023 dropped-ply specimen (designated as number 6-2 in the report 
by Curry, Johnson, and Starnes, 1986) resulted in failure initiation by delamination at the intersec- 
tion of the free edge and the nodal line in the thin section of the laminate at a load of 10,400 lb 
(N x = 3, 467 lb/in. ). The ultimate load of this specimen was 1 1,219 lbs. Specimens with stiffer 
dropped-ply sublaminates delaminated at the drop-off site, and many of the [90 n ], n= 2, 4, 8, 
dropped-ply specimens failed in the end grips. Thus, the failure initiation site for this particular 
[ 9 O 2 ] dropped-ply specimen is somewhat unique. 

A schematic of the finite element modeling scheme for specimen 6-2 is shown in Fig. 2. Solid 
and transition elements are used in the designated shaded region of this figure, which is adjacent 
to the free edge at the nodal line location in the thin section of the laminate. Since symmetry of 
the model about the line y = b/2 is assumed, only one free edge region of solid elements is 
included in the model. The region of solid elements constitutes 5.3% of the planform area in Fig. 
2, and the remaining area is modeled with shell/plate elements. The elements used in the COMET 
software system are: a nine-noded assumed natural strain shell/plate element EX97, a twenty- 
noded serendipity brick element BR20, and a twelve-noded (TR12) and a fifteen-noded (TR15) 
transition element for the connection of the brick elements to the shell/plate elements. 



Fig. 2 Plan view of the [90 2 ] dropped-ply sublaminate specimen 6-2 with approximate locations of 
the nodal lines for the out-of-plane displacement and the region of solid elements indicated. Shell/plate 
elements are used for the remainder of the laminate model. 


Six BR20 elements span the thickness in the solid element region as shown in Fig. 3. Since 
there are sixteen plies in this region, some of the elements span a sublaminate. Effective three- 
dimensional moduli are used to compute element stiffnesses in these cases where one finite ele- 
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ment includes more than one lamina. These effective moduli are calculated by the LAU processor 
in COMET. A theory for three-dimensional modeling is presented by Sun and Li (1988). The 
material properties used for AS4/3502 are given in Table 1. The ply thickness used for the thin 
section is 0.00543 inches and for the thick section the ply thickness is 0.00544 inches. These 
thickness values were determined from measurements of specimen 6-2. 




interfaces between elements 
spanning sublaminates 


5 

4 

3 

2 
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Fig. 3 Six brick elements span the thickness of the sixteen ply laminate in the solid element model- 
ing region. An element occupies a lamina or sublaminate as indicated by the six shaded regions. 


TABLE 1 . Material properties and strength allowables assumed for AS4/3502 
graphite/epoxy unidirectional tape and neat resin 




AS4/3502 graphite/epoxy 

Ei 

126 GPa 

18.5 Msi 

Young’s modulus in fiber direction 

Ez 

11.3 GPa 

1.64 Msi 

Young’s modulus in matrix direction 

e 3 

11.3 GPa 

1.64 Msi 

Young’s modulus in thickness direction 


0.3 


Major Poisson’s ratio in 1-2 plane 


0.3 


Major Poisson’s ratio in 1*3 plane 

^23 

0.3 


Poisson’ ratio in 2-3 plane 

G12 

6.0 Gpa 

0.87 Msi 

Shear modulus in 1-2 plane 

G,3 

6.0 GPa 

0,87 Msi 

Shear modulus in 1-3 plane 

g 23 

3.38 GPa 

0.49 Msi 

Shear modulus in 2-3 plane 



Neat Resin 

E 

3.45 GPa 

0.50 Msi 

Young’s modulus 

u 

0.41 


Poisson’s ratio 

G 

1.34 GPa 

0.195 Msi 

Shear modulus 

Z S1 


Interlaminar Strengths 

93.1 MPa 

13.5 Ksi 

Longitudinal shear strength 

Z S2 

93.1 MPa 

13.5 Ksi 

Transverse shear strength 

Z T 

52.0 MPa 

7.54 Ksi 

Tensile strength 


The load is a uniform axial compressive displacement prescribed at the thick end, and the thin 
end is fixed. Knife edges ran the laminate’s length 0.125 inches inboard of the free edge, and this 
condition is modeled by setting the out-of-plane displacement to zero for the shell elements along 
the knife edge line (see Fig 2). For each of the short ends, the out-of-plane displacement is pre- 
scribed to vanish 0.375 inches into the specimen to simulate the clamped end fixtures of the test. 
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Numerical results are presented for a finite element model consisting of 72 EX97 shell elements, 
180 BR20 solid elements, 96 TR15 transition elements, and 12 TR12 transition elements. The 
total number of active degrees of freedom is 5,263. A total Lagrangian formulation is employed 
along with Newton’s method to solve for the geometrically nonlinear response. 

Interlaminar stresses from the nonlinear finite element analysis are used to compute a dimen- 
sionless delamination index F. This index is defined by 


F = 


fi \ 
' zx ' 


z 51 




x >0 


zz 


( 1 ) 


in which x zz is the interlaminar normal stress component, and x zx and x zy are interlaminar shear 

stress components. The strength allowables Z S1 , Z S2 , and Z T assumed for AS4/3502 are given in 
Table 1. If the index F equals or exceeds unity, then failure by delamination is assumed to have 
initiated. However, at stress risers like the free edge of a laminate where interlaminar stress singu- 
larities are likely to occur in an elasticity solution to effective modulus models, the peak value of 
Fis meaningless. For finite element approximations to singular elasticity solutions, peak values of 
F at the stress riser increase monotonically with increasing mesh refinement. Hence the applica- 
tion of a point stress criterion like Eq. (1) has to be modified when singular behavior is expected. 
In the case of singular behavior, interlaminar stresses from the finite element analysis may be 
average over some distance from the stress riser and these average values used in Eq. (1). The 
application of Eq. (1) with averaged stresses was used by Brewer and Lagace (1988) as a delami- 
nation criterion. We (Ddvila and Johnson, 1992) used Eq. (1) with stresses evaluated at a distance 
d from the ply drop-off as a delamination criterion. Singular behavior is expected in the analysis 
of specimen 6-2, so the distributions of F from the finite element analysis can be expected to 
exceed unity at loads near the experimental failure initiation load. 

The results from the finite element analysis at 3,433 lb/in. compression are shown in Figs 4-8. 
This is the load increment from the analysis just below the experimental failure load of 3,467 lbs/ 
in. The delamination index is plotted over the interface between the -45 degree ply and the +45 
degree ply in Fig. 4, which corresponds to interface 1 in Fig. 3. This interface had the largest mag- 
nitude of F of the six interfaces where F was evaluated. For the distribution shown in Fig. 4, F max 

= 1.073 at y = 0 (free edge) and x = 6.778 inches. Also, F = 1 atx = 6.778 when y;= 0.0063 inches 
in Fig. 4. At a load of 3,516 lbs/in, which is the load increment in the analysis just above the 
experimental failure load, F= 1 at y = 0.0095 inches. Thus, by interpolation the characteristic dis- 
tance d inward from the free edge for which Eq. (1) would be used to predict delamination initia- 
tion is about 0.0076 inches, or d is equal to approximately 1.4 ply thicknesses. 

The individual interlaminar stress components which determine the F-distribution in Fig. 4 
are shown in Figs. 5-7. Examination of these figures show that the stress component t zx , the inter- 
laminar shear stress component tangent to the free edge, is the largest contributor to F. Since x xz is 
substantially larger than x zz and x zy at the free edge in this interface (compare Fig. 6 to Figs. 5 and 
7), it is likely that t h initiates delamination in the -45/+45 interface at the intersection of the 
nodal line with the free edge in the thin section of the dropped-ply laminate. 

Finally, the influence of the geometric nonlinearity on the interlaminar stress distribution can 
be illustrated by computing the delamination index from a linear analysis at the same free edge 
location as used for the geometrically nonlinear analysis for specimen 6-2. The linear analysis 
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results differs from the nonlinear results in that only one nodal line in the out-of-plane displace- 
ment occurs rather than the two nodal lines in the nonlinear response. The one nodal line in the 
linear analysis occurs near the ply drop-off, and no nodal line occurs in the thin section. The dis- 
tribution of the delamination index F at interface 2 (see Fig. 3) from the linear analysis is shown 
in Fig. 8 at a compressive load of 3,433 lb/in. Interface 2 had a larger peak value of F than inter- 
face 1. Comparing Figs. 4and 8, it can be seen that a linear analysis substantially underpredicts 
the interlaminar stress intensity for specimen 6-2. 


w 



Fig. 4 Distribution of the delamination index F in the -45/+45 interfacial plane (interface 1 in Fig. 
3) of specimen 6-2 from the nonlinear analysis at a compressive load of 3,433 lb/in. Free edge at y = 
0 and knife edge support at y = 0.125 in. 


Graduate research assistant 

Ms. F. Aylin Barlas is the graduate assistant on this project, and this work is the basis for her 
master’s thesis in Aerospace Engineering. 
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Fig. 5 Distribution of the interlaminar normal stress component in the -45/+4S interfacial plane 
(interface 1 in Fig. 3) of specimen 6-2 from nonlinear analysis at a compressive load compressive 
load of 3,433 lb/in. Free edge at y = 0 and knife edge support at y = 0. 125 in. 



Fig. 6 Distribution of the interlaminar shear stress component tangent to the free edge in the -45/+45 
interfacial plane (interface 1 in Fig. 3) of specimen 6-2 from nonlinear analysis at a compressive load of 
3,433 lb/in. Free edge at y = 0 and knife edge support at y = 0.125 in. 
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Fig. 7 Distribution of the interlaminar shear stress component normal to the free edge in the -45/ 
+45 interfacial plane of specimen 6-2 from nonlinear analysis at a compressive load of 3,433 lb/in. 
Free edge at y = 0 and knife edge support at y = 0.125 in. 



Fig. 8 Distribution of the delamination index F in the +45A45 interfacial plane (interface 2 in Fig. 3) 
of specimen 6-2 from linear analysis at a compressive load of 3,433 Ib/in. Free edge at y = 0 and knife 
edge support at y = 0.125 in. 
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Stiffener Crippling Initiated by Delamination 


Background and objective 

Graphite-epoxy stiffeners tested in compression (Bonanni, Johnson, and Starnes, 1991) 
showed evidence of failure initiation in postbuckling by delamination at the free edges of the 
flanges. This mode of failure initiation occurred in I- and 7-section specimens with flange width to 
thickness ratios in the range from ten to twenty. The analyses conducted in the Bonanni study 
modeled the stiffeners as a branched shell using the STAGS computer code. Since shell theories 
are based on plane stress and neglect through-the-thickness stress components, failure initiated by 
interlaminar stresses at the free edges of the flanges due to bending in postbuckling could not be 
studied with the STAGS model. 

The objective of this project is to develop a computational model of the stiffener specimens 
that includes the capability to predict the interlaminar stress response at the flange free edge in 
postbuckling. 

I-section specimen 

Computational models were developed for the I-section specimen designated 1-2 by Bonanni, 
et al. (1991). This specimen had an axial gage length of 5.944 inches with an additional inch at 
each end of the specimen potted into fixtures. The cross section and measured dimensions are 
shown in Fig. 9 below. Specimens were fabricated from AS4/3502 graphite/epoxy unidirectional 
tape with material properties listed in Table 1. The wall construction was eight plies of tape with a 
[±45/0/90] s stacking sequence in the web and two diagonally opposite flanges. The stacking 

sequence in the two other flanges is [±45/0/90 2 /0/±45] because of the fabrication process 

with tape sublaminates. Crippling for specimen 1-2 initiated due to a local material compressive 
strength failure in the comer where the flanges and web meet and not due to flange free edge 
delamination. However, the detailed STAGS results reported by Bonanni et al. (1991) for speci- 
men 1-2 provides a good benchmark to assess current modeling techniques. Once current model- 
ing techniques are established for this benchmark case, then this experience will guide the 
modeling of specimens that did initiate crippling by flange free edge delamination. 
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Dimensional data 
t = 0.0455 inches 

bf= 1.262 inches. 
b w = 1 .247 inches 
r= 0.125 inches 


Computational modeling 

A geometrically nonlinear finite element model is being developed using the COMET soft- 
ware system. The model consists of nine-noded EX97 shell/plate elements, transition elements 
TR12 and TR15, and twenty-noded solid elements BR20. The solid elements are used in a local 
region near the flange free edge to estimate the interlaminar stresses. They are connected to transi- 
tion elements which in tum are connected to the shell elements. Shell/plate elements are used in 
the majority of the stiffener. Thus, the solid elements are only used in regions where three-dimen- 
sional effects are important, and this is a significant modeling strategy because a stiffener modeled 
with all solid elements is not warranted based on the mechanics and is prohibitively large compu- 
tationally. This issue of very large degree-of-freedom numerical models is compounded in nonlin- 
ear analysis because of the iterative nature of the solution procedure. 

Two shell models were examined initially; i.e., no solid element subdomains were included. 
Both models used EX97 plate elements in the four flanges and web. It is at the junctions of the 
web to flanges where the two models differed. As shown in Fig. 9, the actual specimen had zero- 
degree tape inserted into the junctions to fill the void created by sublaminate fabrication the spec- 
imens. This creates a stiff insert in the comer junctions. Model I use E210 beam elements to rep- 
resent the stiff insert, and Model II used narrow-width plate elements with a zero-degree 
unidirectional layup properties to represent the stiff insert. A bifurcation buckling analysis was 
performed for the purely compressive prebuckling equilibrium configuration. Model I had 3,925 
active degrees of freedom and Model II had 4,515 active degrees of freedom. The results for the 
load-end shortening stiffness in prebuckling equilibrium K pb and the bifurcation buckling load P cr 
are shown in Table 2 along with experimental results. The results attained with Model II are 
slightly better than with Model I. However, Model II had longer CPU times to converge relative to 
Model I because of its increased number of degrees of freedom relative to Model I. Model I was 
close enough to the experiment to select it for further nonlinear analysis. 
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Table 2: Prebuckling stiffness AT pb and bifurcation buckling load P cr from two finite element 

models of specimen 1-2 compared to experiment 


Xp b , lb/in. P QV lb. 


Model I a 

282,800 

2,931 

Model II b 

302,193 

2,708 

Experiment (Bonanni, et al, 1991) 

296,450 

2,786 


a. For shell elements E ( = 17 Msi, and for beam elements E = 10 Msi. 

b. For all shell elements Ej = 18.5 Msi. 


Models I and II were used in a geometrically nonlinear analysis by including an initial shape 
imperfection in the pattern of the buckling mode with a maximum amplitudes of one-tenth of the 
flange thickness. The load-end shortening plot from these analyses is shown in Fig. 10 along with 
the experimental results. Both models agree very well with the experimental results. Model II 
results are only shown to a compressive load magnitude slightly greater than 6,000 lbs, since it 
converges more slowly than Model I. These nonlinear results again indicate that the shell/plate 
Model I is adequate to capture the load-shortening response of specimen 1-2. 

Adding a reasonable number of solid elements in a subdomain at the free edge of a flange for 
Model I increased the degrees of freedom from 3,925 in the shell/beam model to over 5,500 for 
the solid/shell/beam model. Since we are only interested in computing interlaminar stresses near 
the experimental failure load, and the shell/beam model is apparently sufficient to capture the geo- 
metrically nonlinear response, we propose to do variable complexity modeling. That is, a nonlin- 
ear analysis is carried out on the shell/beam model (3,925 degrees of freedom) to a load near 
failure initiation in the experiment. Then, the complexity of the model is changed to include the 
subdomain of solid elements, an additional 1,500 degrees of freedom or so. Since Newton’s 
method is based on a good initial guess to the solution, shell kinematics may be employed to esti- 
mate the three displacements of each node in the solid element subdomain. With this initial guess, 
Newton’s method is initiated for the complex solid/shell/beam model. Hopefully, we can converge 
at the last successful load step of the shell/beam model with only one sequence of Newton’s itera- 
tions for the complex model at this same load. This strategy requires more programming in the 
COMET software system than attempted to date. Displacement data files have to created from the 
shell/beam model that are applicable to the solid/shell/beam model, so that a restart of the 
NL_STATIC_1 procedure on the complex model can be initiated. We have worked-out the pro- 
gramming details on simpler problems involving variable Complexity modeling, and are in the 
process of applying these programming steps to the I-section model. 
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failure 



Axial shortening in inches 


Fig. 10 Load-end shortening for specimen 1-2 compared to two finite element models of its response. 


Graduate research assistant 

Ms. F. Aylin Barlas is the graduate assistant on this project, and this work is the basis for her 
master’s thesis in Aerospace Engineering. 
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Pressure Pillowing of an Orthogonally Stiffened Cylindrical Shell 
Objective and background 

This project is motivated by the application of advanced composite materials in the construc- 
tion of fuselage primary structure for transport aircraft. We consider a laminated composite mate- 
rial cylindrical shell internally stiffened by a regular arrangement of stringers and frames, or rings, 
and subjected to internal pressure and axial load. The objective is to determine the distribution of 
the interacting loads between the stiffeners and the shell wall; particularly the load transfer at the 
stiffener crossing point. 

Summary of approach 

For the mathematically infinite length circular cylinder, a structural repeating unit is defined 
that consists of a portion of the shell wall centered over the stringer and frame intersection. See 
Fig. 11. The shell is modeled with Sander’s theory, and the stiffeners are modeled with Vlasov’s 
thin-walled bar theory extended to composite material construction. The principle of virtual work 
is employed for each structural element, and the assumed 3X1 displacement vector functions n 
for each element is the sum of a principal part n p , a particular part Up and a Fourier Series part 
»f S . In the shell, for example, the principal part of the displacement vector contains exponentially 

decaying functions which account for the jump in the stress and moment resultants at the stiffener 
connections, the particular part accounts for stretching under axial load, and the Fourier Series 
part is periodic in the frame and stringer spacing. Inter-element continuity is enforced by aug- 
menting the virtual work functional with integrals of Lagrange multiplier functions times varia- 
tions in the displacement constraints. The Lagrange multipliers are interpreted as the components 
of the line loads between the stiffeners and shell, and they are expanded in a Fourier Series as 
well. Augmented virtual work leads to a set of equilibrium equations for the displacement and 
Lagrange multiplier coefficients. The displacement constraint equations are enforced by the van- 
ishing of the inner product of these equations with the variations in the Lagrange multiplier func- 
tions. The complete set of equations for the displacement and Lagrange multiplier coefficients is 
indeterminate, but can be solved by first eliminating the displacement sub-vector and then solving 
for the Lagrange multiplier sub-vector. 



Fig. 1 1 Repeating unit of an orthogonally stiffened cylindrical shell subject to internal pressure p and 
axial load P. 
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Preliminary results 

We are comparing the results of our model to those published by Wang and Hsu (1985). Wang 
and Hsu give results for a balanced, symmetrically laminated skin with a 1 17.5 inch radius, a 20 
inch frame spacing (L = 20), a 5.8 inch stringer spacing (b = 5.8), and for an internal pressure p of 
10 psi. They include closed-end pressure vessel effect which implies that the axial load P = 
3,407.5 lb on a generic repeating unit. The results presented below include the Fourier series por- 
tion of the assumed displacement fields and the particular part to account for axial stretching. The 
principal part of the displacement field that accounts for jump conditions in the shell in-plane 
shear and Kirchhoff shear at the stiffener attachment lines are not included. The Fourier Series 
expansions were truncated at twenty-four harmonic terms in the axial and circumferential direc- 
tions. 

The distributions of the out-of-plane displacements of the skin from our solution are shown in 
Figs. 12 and 13. These displacements compare very well to those presented by Wang and Hsu (see 
their Fig. 5). Axial strain distributions on the inner and outer surfaces of the skin at an internal 
pressure of 10 psi are shown in Fig. 14, and the circumferential strain distributions on the surfaces 
of the skin are shown in Fig. 15. These strain distributions do not agree with those presented by 
Wang and Hsu (see their Figs. 6 and 7). For example, the axial strain distribution in the circumfer- 
ential direction from our results (Fig. 14) starts at about 2000)18 at x = 0 and y =0 and remains the 
same value to x = 0 and y = b/2. The corresponding strain distribution from Wang and Hsu’s 
results is about 2700pe at x =y = 0 and decreases to about 45|xe at x = 0 and y = b/2. The decrease 
in axial strain from 2700)18 to 45)18 is attributed to the present of the stringer at y = b/2, and it is 
clear our results do not capture this influence. We believe that including the principal part of the 
shell displacement field will help alleviate this discrepancy, and we are also investigating if the 
particular part of the displacement field that accounts for stretching due to axial load is suffi- 
ciently general. The trends in the remaining strain distributions shown in Figs. 14 and 15 are sim- 
ilar to Wang and Hsu’s results but do not agree quantitatively. 

An interesting result not published by Wang and Hsu, is the distribution of the interacting line 
loads between the skin and stiffeners. The normal and tangential components of the line load 
between the stringer and skin are shown in Figs. 16 and 17, respectively. The corresponding com- 
ponents between the frame and skin are shown in Figs. 18 and 19. The normal components of 
these line loads exhibit delta function behavior, and the tangential components exhibit step func- 
tion behavior, at the stiffener crossing point. Although the results are not presented here, it was 
observed that increasing the number of terms in the Fourier series, within a range limited by com- 
puter resources, caused increasing peak magnitudes of the normal line loads at the crossing point. 
This interesting result suggests that we should include displacement fields to account for singular- 
ities in the force and moment resultants in the skin, and point forces at the crossing point in the 
Lagrange multiplier distribution. 

Graduate research assistant 

Naveen Rastogi is the graduate assistant on this project, and this work is the basis of his dis- 
sertation for the Doctor of Philosophy degree in Aerospace Engineering. 
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Fig. 12 Axial distributions of the out-of-plane displacement of the skin at three circumferential positions 
at an internal pressure of 10 psi. 
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Fig. 13 Circumferential distributions of the out-of-plane displacement of the skin at three axial loca- 
tions at an internal pressure of 10 psi. 
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Fig. 14 Axial normal strain distributions on the skin surfaces at an internal pressure of 10 psi. 
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Fig. 15 Circumferential strain distributions on the skin surfaces at an internal pressure of 10 psi. 
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Fig. 16 Normal load intensity between the stringer and skin at an internal pressure of 10 psi. 
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Fig. 17 Tangential load intensity between stringer and skin at an internal pressure of 10 psi. 
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Fig. 18 Normal load intensity between frame and skin at an internal pressure of 10 psi. 
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Fig. 19 Tangential load intesity between frame and skin at an internal pressure of 10 psi. 
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